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In the bulk state, the Ising FCC antiferromagnet is fuhy frustrated and is known 
to have a very strong first-order transition. In this paper, we study the nature of this 
phase transition in the case of a thin film, as a function of the film thickness. Using 
Monte Carlo (MC) simulations, we show that the transition remains hrst order down 
to a thickness of four FCC cells. It becomes clearly second order at a thickness of 
two FCC cells, i.e. four atomic layers. It is also interesting to note that the presence 
of the surface reduces the ground state (GS) degeneracy found in the bulk. For 
the two-cell thickness, the surface magnetization is larger than the interior one. It 
undergoes a second-order phase transition at a temperature Tc while interior spins 
become disordered at a lower temperature T^. This loss of order is characterized by 
a peak of the interior spins susceptibility and a peak of the specific heat which do 
not depend on the lattice size suggesting that either it is not a real transition or it 
is a Kosterlitz-Thouless nature. The surface transition, on the other hand, is shown 
to be of second order with critical exponents deviated from those of pure 2D Ising 
universality class. We also show results obtained from the Green's function method. 
Discussion is given. 

PACS numbers: 75.10.-b General theory and models of magnetic ordering ; 75.40.Mg Nu- 
merical simulation studies ; 75.70.Rf Surface magnetism 
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I. INTRODUCTION 



This paper deals with the question whether or not the phase transition known in the bulk 
state changes its nature when the system is made as a thin film. In a recent work, we have 
considered the case of a bulk second-order transition. We have shown that under a thin film 
shape, i.e. with a finite thickness, the transition shows effective critical exponents whose 
values are between 2D and 3D universality classes.- If we scale these values with a function 
of thickness as suggested by Fisher- we should find, as long as the thickness is finite, the 2D 
universality class. 

In this paper, we study the effect of the film thickness in the case of a bulk first-order 
transition. The question to which we would like to answer is whether or not the first order 
becomes a second order when reducing the thickness. For that purpose we consider the 
face-centered cubic (FCC) Ising antiferromagnet. This system is fully frustrated with a very 
strong first-order transition. 

On the one hand, effects of the frustration in spin systems have been extensively in- 
vestigated during the last 30 years. In particular, by exact solutions, we have shown that 
frustrated spin systems have rich and interesting properties such as successive phase transi- 
tions with complicated nature, partial disorder, reentrance and disorder lines.-i^ Frustrated 
systems still challenge theoretical and experimental methods. For recent reviews, the reader 
is referred to Ref. |5. 

On the other hand, physics of surfaces and objects of nanometric size have also attracted 
an immense interest. This is due to important applications in industry.-'^'^ In this field, 
research results are often immediately used for industrial applications, without waiting for a 
full theoretical understanding. An example is the so-called giant magneto-resistance (GMR) 
used in data storage devices, magnetic sensors, ..M-^^^ In parallel to these experimental 
developments, much theoretical effort has also been devoted to the search of physical mech- 
anisms lying behind new properties found in nanometric objects such as ultrathin films, 
ultrafine particles, quantum dots, spintronic devices etc. This effort aimed not only at 
providing explanations for experimental observations but also at predicting new effects for 
future experiments.—!^ 

The above-mentioned aim of this paper is thus to investigate the combined effects of 
frustration and film thickness which are expected to be interesting because of the symmetry 



reduction. As said above, the bulk FCC Ising antiferromagnet is fully frustrated because it 
is composed of tetraliedra whose faces are equilateral triangles. The antiferromagnetic (AF) 
interaction on such triangles causes a full frustration.^ The bulk properties of this material 
have been largely studied as we will show below. In this paper, we shall use the recent high 
precision technique called "Wang-Landau" flat histogram Monte Carlo (MC) simulations 
to identify the order of the transition. We also use the Green's function (GF) method for 
qualitative comparison. 

The paper is organized as follows. Section II is devoted to the description of the model. 
We recall there properties of the 3D counterpart model in order to better appreciate prop- 
erties of thin films obtained in this paper. In section III, we show our results obtained by 
MC simulations on the order of the transition. A detailed discussion on the nature of the 
phase transition is given. In the regime of second-order transition, we show in this section 
the results on the critical exponents obtained by MC flat histogram technique. Section IV is 
devoted to a study of the quantum version of the same model by the use of the GF method. 
Concluding remarks are given in section V. 



II. MODEL AND GROUND STATE ANALYSIS 



It is known that the AF interaction between nearest-neighbor (NN) spins on the FCC 
lattice causes a very strong frustration. This is due to the fact that the FCC lattice is 
composed of tetrahedra each of which has four equilateral triangles. It is well-known^ that it 
is impossible to fully satisfy simultaneously the three AF bond interactions on each triangle. 
In the case of Ising model, the GS is infinitely degenerate for an infinite system size: on each 
tetrahedron two spins are up and the other two are down. The FCC system is composed 
of edge-sharing tetrahedra. Therefore, there is an infinite number of ways to construct the 
infinite crystal. The minimum number of ways of such a construction is a stacking, in one 
direction, of uncorrelated AF planes. The minimum GS degeneracy of a Li^ FCC-cell system 
(L being the number of cells in each directions), is therefore equal to 3 x 2^^ where the 
factor 3 is the number of choices of the stacking direction, 2 the degeneracy of the AF spin 
configuration of each plane and 2L the number of atomic planes in one direction of the FCC 
crystal (the total number of spins is = 4L^). The GS degeneracy is therefore of the order 
of 2^^'^ . Note that at finite temperature, due to the so-called "order by disorder",—*^ the 



spins will choose a long-range ordering. In the case of AF FCC Ising crystal, this ordering is 
an alternate stacking of up-spin planes and down-spin planes in one of the three direction. 
This has been observed also in the Heisenberg well as in other frustrated systems.— 

The phase transition of the bulk frustrated FCC Ising antiferromagnet has been found 
to be of the first order.— '^'^'^'^ Note that for the Heisenberg model, the transition is also 
found to be of the first order as in the Ising case.— Other similar frustrated antiferromag- 
nets such as the HCP XY and Heisenberg antiferromagnets^^ and stacked triangular XY 
and Heisenberg antiferromagnets2^>2I show the same behavior. 

Let us consider a film of FCC lattice structure with (001) surfaces. The Hamiltonian is 
given by 

'H = -Y^ Jijai ■ Gj (1) 

where Oi is the Ising spin at the lattice site i, indicates the sum over the NN spin pairs 

Gi and Oy 

In the following, the interaction between two NN on the surface is supposed to be AF 
and equal to J^. All other interactions are equal to J = —1 for simplicity. Note that in a 
previous paper,— we have studied the case of the Heisenberg model on the same FCC AF 
film as a function of Jg- 

For Ising spins, the GS configuration can be determined in a simple way as follows: we 
calculate the energy of the surface spin in the two configurations shown in Fig. [1] where the 
film surface contains spins 1 and 2 while the beneath layer spins 3 and 4. In the ordering 
of type I [Fig. [11(a)], the spins on the surface [xy plane) are antiparallel and in the ordering 
of type II [Fig. [11(b)] they are parallel. Of course, apart from the overall inversion, for 
type I there is a degenerate configuration by exchanging the spins 3 and 4. To see which 
configuration is stable, we write the energy of a surface spin for these two configurations 

Ei = -A\Js\ 

= 4|J,|-4|J| (2) 

One sees that < Ejj when Js/J > 0.5. In the following, we study the case Jg = J = —I 
so that the GS configuration is of type I. 



(a) 



(b) 



FIG. 1: The ground state spin configuration of tlie FCC cell at the film surface: a) ordering of 
type I for Jg < — 0.5| J| b) ordering of type II for Jg > — 0.5| J|. 



In this paragraph, we show the results obtained by MC simulations with the Hamiltonian 
([I]) using the high-precision Wang-Landau flat histogram technique.— Wang and Landau 
recently proposed a MC algorithm for classical statistical models. The algorithm uses a 
random walk in energy space in order to obtained an accurate estimate for the density 
of states p{E) which is defined as the number of spin configurations for any given E. This 
method is based on the fact that a fiat energy histogram H{E) is produced if the probability 
for the transition to a state of energy E is proportional to p{E)^^. At the beginning of the 
simulation, the density of states (DOS) is set equal to one for all possible energies, p{E) = 1. 
We begin a random walk in energy space (E) by choosing a site randomly and flipping its 
spin with a probability proportional to the inverse of the momentary density of states. 
In general, if E and E' are the energies before and after a spin is flipped, the transition 
probability from E to E' is 



Each time an energy level E is visited, the DOS is modified by a modification factor / > 
whether the spin flipped or not, i.e. p{E) p{E)f. At the beginning of the random walk, 
the modification factor / can be as large as ~ 2.7182818. A histo gram H{E) records 
how often a state of energy E is visited. Each time the energy histogram satisfies a certain 
"flatness" criterion, / is reduced according to / — v7 H{E) is reset to zero for all 
energies. The reduction process of the modification factor / is repeated several times until 
a final value /gnai which close enough to one. The histogram is considered as flat if 



III. MONTE CARLO RESULTS 



p{E E') = min [p{E)/p{E'), 1] . 



(3) 



H{E) > x%.{H{E)) 



(4) 



for all energies, where x% is chosen between 70% and 95% and {H{E)) is the average 
histogram. 

The thermodynamic quantities'"^ can be evaluated by 

(E") = WEy{E)exp{~E/kBT) 

^ E 

_ {E^) - {Ey 

(M") = WMy{E)expi-E/kBT) 

^ E 

^ ~ w 

where Z is the partition function defined by 

Z = Y^p{E)^xp{-ElkBT) (5) 

E 

The canonical distribution at any temperature can be calculated simply by 

P(E, T) = ^p{E) eM-E/kBT) (6) 

In this work, we consider a energy range of interest^*^ (i?min, -E'max)- We divide this 
energy range to R subintervals, the minimum energy of each subinterval is -Emin ^r i = 
1, 2, R, and maximum of the subinterval i is -E^ax — -^min + ^A-^', where AE can be chosen 
large enough for a smooth boundary between two subintervals. The Wang-Landau algorithm 
is used to calculate the relative DOS of each subinterval {E^^^^, -E^ax) with the modification 
factor /final = exp(10~^) and flatness criterion x% = 95%. We reject the suggested spin flip 
and do not update p{E) and the energy histogram H{E) of the current energy level E if the 
spin-flip trial would result in an energy outside the energy segment. The DOS of the whole 
range is obtained by joining the DOS of each subinterval (-Emin + ^-^5 -^max ~ AE). 

The film size used in our present work is L x L x where L is the number of cells 
in X and y directions, while is that along the z direction (film thickness). We use here 
L = 30,40, 150 and A^^ = 2, 4, 8, 12. Periodic boundary conditions are used in the xy 
planes. Our computer program was parallelized and run on a rack of several dozens of 64-bit 
CPU. I J| = 1 is taken as unit of energy in the following. 

Before showing the results let us adopt the following notations. Sublattices 1 and 2 of 
the first FCC cell belongs to the surface layer, while sublattices 3 and 4 of the first cell 



FIG. 2: Bulk energy vs T for L = iV^ = 12. 
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FIG. 3: Bulk energy histogram for L = Nz = 12 with periodic boundary conditions in all three 
directions (a) and without PBC in z direction (b). The histogram was taken at the transition 
temperature Tc indicated on the figure. 

belongs to the second layer [see Fig. [Uj^a)]. In our simulations, we used A^^ FCC cells, i.e. 
2Nz atomic layers. We used the symmetry of the two film surfaces. 

A. Crossover of the phase transition 

As said earlier, the bulk FCC antiferromagnet with Ising spins shows a very strong first- 
order transition. This is seen in MC simulation even with a small lattice size as shown in 
Fig. El 

Our purpose here is to see whether the transition becomes second order when we decrease 
the film thickness. As it turns out, the transition remains of first order down to A^^ = 4 as 
seen by the double-peak energy histogram displayed in Fig. HI Note that we do not need to 
go to larger L, the transition is clearly of first order already at L = 40. 

In Fig. Owe plot the latent heat Ai? as a function of thickness A^^. Data points are well 
fitted with the following formula 




FIG. 4: Energy histogram for L = 20, 30, 40 with film thickness A^^ = 4 (8 atomic layers) at 
T = 1.8218, 1.8223, 1.8227, respectively. 
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FIG. 5: The latent heat IS.E as a function of thickness A^^ 



AE = A- ^ 



where d = 3 is the dimension, A = 0.3370, B = 3.7068, C = -0.8817. Note that the term 
N^~^ corresponds to the surface separating two domains of ordered and disordered phases 
at the transition. The second term in the brackets corresponds to a size correction. As seen 
in Fig. [5l the latent heat vanishes at a thickness between 2 and 3. This is verified by our 
simulations for Nz = 2. For Nz = 2 we find a transition with all second-order features: no 
discontinuity in energy (no double-peak structure) even when we go up to L = 150. 

Before showing in the following the results of A^^ = 2, let us discuss on the crossover. In 
the case of a film with finite thickness studied here, it appears that the first-order character 
is lost for very small A'^^. A possible cause for the loss of the first-order transition is from the 
role of the correlation in the film. If a transition is of first order in 3D, i. e. the correlation 
length is finite at the transition temperature, then in thin films the thickness effect may 
be important: if the thickness is larger than the correlation length at the transition, than 
the first-order transition should remain. On the other hand, if the thickness is smaller than 



FIG. 6: Energy versus temperature T for L = 120 with film thickness A'^^ = 2. 
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FIG. 7: Sublattice magnetization for L = 120 with film thickness iV^ = 2. 

that correlation length, the spins then feel an "infinite" correlation length across the film 
thickness resulting in a second-order transition. 

B. Film with 4 atomic layers (A^^ = 2) 

Let us show in Fig. [6] and Fig. [7] the energy and the magnetizations of sublattices 1 and 
3 of the first two cells with L = 120 and = 2. 

It is interesting to note that the surface layer has larger magnetization than that of the 
second layer. This is not the case for non frustrated films where the surface magnetization 
is always smaller than the interior ones because of the effects of low-lying energy surface- 
localized magnon modes.— One explanation can be advanced: due to the lack of neighbors 
surface spins are less frustrated than the interior spins. As a consequence, the surface spins 
maintain their ordering up to a higher temperature. 

Let us discuss on finite-size effects in the transitions observed in Figs. [H] and This is 
an important question because it is known that some apparent transitions are artifacts of 
small system sizes. 
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FIG. 8: Specific heat are shown for various sizes L as a function of temperature. 
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FIG. 9: Susceptibihties of sublattice 1 (a) and 3 (b) are shown for various sizes L as a function of 
temperature. 



To confirm further the observed transitions, we have made a study of finite-size effects 
on the layer susceptibihties by using the Wang-Landau technique described above.— 

We observe that there are two peaks in the specific heat: The first peak at Ti ~ 1.927, 
corresponding to the vanishing of the sublattice magnetization 3, does not depend on the 
lattice size while the second peak at T2 ^ 1.959, corresponding to the vanishing of the 



sublattice magnetization 1, does depend on L. Both histograms taken at these temperatures 
and the near-by ones show a gaussian form indicating a non first-order transition [see Fig. 
HO]. 
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FIG. 10: Energy histograms for L = 120 with film thickness Nz = 2 at the two temperatures 
(indicated on the figure) corresponding to the peaks observed in the specific heat. See text for 
comment. 

The fact that the peak at Ti does not depend on L suggests two scenarios: 

i) the peak does not correspond to a real transition, since there exist systems where Cy 
shows a peak but we know that there is no transition just as in the case of ID Ising chain, 

ii) the peak corresponds to a Kosterhtz-Thouless transition. To confirm this we need to 
check carefully several points such as the behavior of the correlation length etc. This is a 
formidable task which is not the scope of this work. 

Whatever the scenario for the origin of the peak at Ti, we know that the interior layers 
are disordered between Ti and T2, while the two surface layers are still ordered. Thus, 
the transition of the surface layers occurs while the disordered interior spins act on them 
as dynamical random fields. Unlike the true 2D random-field Ising model which does not 
allow a transition at finite temperature,— the random fields acting on the surface layer are 
correlated. This explains why one has a finite-T transition here. Note that this situation 
is known in some exactly solved models where partial disorder coexists with order at finite 
However, it is not obvious to foresee what is the universality class of the transition 
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FIG. 11: The maximum value of Vi = {E) — {ME) / (M) versus L in the In — In scale. The slope 
of this straight line gives 

at T2. The theoretical argument of Capehart and Fisher- does not apply in the present 
situation because one does not have a single transition here, unlike the case of simple cubic 
ferromagnetic films studied before.^ So, we wish to calculate the critical exponents associated 
with the transition at T2. 

The exponent u can be obtained as follows. We calculate as a function of T the magneti- 
zation derivative with respect to p = {kBT)-^: Vi = ((InM)') = {E) - {ME) / (M) where E 
is the system energy and M the sublattice order parameter. We identify the maximum of Vi 
for each size L. From the finite-size scaling we know that V^"^^ is proportional to L^/'^.^^ We 
show in Fig. [TT]the maximum of Vi versus InL for the first layer. We find u = 0.887± 0.009. 
Now, using the scaling law t^™^^ oc L'^^'^, we plot In;^;™'^^ versus InL in Fig. [T21 The ratio 
of the critical exponents 7/// is obtained by the slope of the straight line connecting the 
data points of each layer. From the value of u we obtain 7 = 1.542 ± 0.005. These values 
do not correspond neither to 2D nor 3D Ising models (72D = 1-75, 1^20 = 1, 73D = 1.241, 
i^3D = 0.63). We note however that, if we think of the weak universality where only ratios 
of critical exponents are concerned,-® then the ratios of these exponents 1/u = 1.128 and 
7/z/ = 1.739 are not far from the 2D ones which are 1 and 1.75, respectively. 



IV. GREEN'S FUNCTION METHOD 



We consider here the same FCC system but with quantum Heisenberg spins. To compare 
the results with the Ising case, we add in the Hamiltonian an Ising-like anisotropy interaction 
term. In addition, this term avoids the absence of long-range order of isotropic non Ising 
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FIG. 12: Maximum sublattice susceptibility t^™^^ versus L in the In — In scale. The slope of this 
straight line gives 

spin model at finite temperature (T) when the film thickness is very small, i.e. quasi two- 
dimensional system.— The Hamiltonian is given by 

^ = - H J^j^i ■ - 51 kjS^S- (8) 

where S, is the Heisenberg spin at the lattice site i, indicates the sum over the NN spin 
pairs Si and S^. Jij and lij are antiferromagnetic (negative). Note that in the laboratory 
coordinates, the up spins have 5*^ > while the down spins have < 0. 

We can rewrite the Hamiltonian ([8]) in the relative local spin coordinates as 

^ = - E Aj\\{cose,,-l){stSl + SiSr) 

<i,j> 

+ ^ (cos% + 1) [S+SJ- + SrS+) 

+ \ sin e,, {St + Si)S^^-\ sin {Sj + Sj) 

+ COS^^fsA - E k,COs0^3StS| (9) 
J <i,i> 

where 9ij is the angle between two NN spins. Note that in the above expression, we have 
transformed all S"^ > 0, the relative spin orientation of each spin pair is now expressed by 
9ij. In a collinear spin configuration such as those shown in Fig. [T], cos 6ij = — 1 and 1 for 
antiparallel and parallel pairs, respectively, while sin 6ij = 0. In non collinear structures, the 
calculation is more complicated. The general GF method for non collinear spin configuration 
has been proposed elsewhere.— In the present study, one has a collinear spin configuration 



shown in Fig. [T] because of the Ising-hke anisotropy. We define two double-time GF by 

Gi,it,t') = (10) 



(11) 



where i and j belong to the up-spin sublattice, m to the down-spin one. In the case of 



thin films, the reader is referred to Refs. U8 



33 



34 for a general formulation. We describe 



here only the main steps: we first write the equations of motion for Gij{t,t') and Fmj{t,t') 
and we next neglect higher-order correlations by using the Tyablikov decoupling scheme^ 
which is known to be valid for exchange terms,— and then we introduce the following Fourier 
transforms in the xy plane 
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(12) 



(13) 



where lo is the spin-wave frequency, h^y denotes the wave-vector parallel to xy planes, Rj is 
the position of the spin at the site i, n and n' are respectively the indices of the layers where 
the sites i (or m) and j belong to. One has n,n' = 1,2, ...,2Nz. The integral over is 
performed in the first Brillouin zone in the xy reciprocal plane whose surface is A. Finally, 
one obtains for all layers the following matrix equation 



M (cu) g = u. 



(14) 



Note that though n runs from 1 to 2Nz, the matrix M has the dimension of AN^ x 4Nz 
because for each n there are two functions g{n,n') and f{n,n'). In the above equation, g 
and u are the column matrices of dimension 4Nz which are defined as follows 
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where for the spin configuration of type I [Fig. [T] (a)] one has 



4± 



in which, Z = 4 is the number of in-plane NN, and 
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4cos — cos — 



(16) 



Bn 


= -2 J„ (S^) (Z7) 


(17) 




= +4J.,_,(5:)cos^ 


(18) 




= -4 J„,„_i (S"^) cos ^ 


(19) 




= +4J„,.+i(5:)cos^ 


(20) 


D- 


= -4 J„,„+i (5^) cos ^ 


(21) 



Here, for compactness we have used the following notations: 

i) Jn and /„ are the in-plane interactions. In the present model Jn is equal to J = — 1. 
All /„ are set to be /(< 0). 

ii) Jn,n±i are the interactions between a spin in the n-th layer and its neighbor in the 
{n ± l)-th layer. Here, we take Jn,n±i = —1- Of course, Jn,n-i = if n = 1, Jn,n+i = if 
n = 2N,. 

Now, solving det(M) = |M| = 0, we obtain the spin-wave spectrum u of the present 
system. For each k^y there are iN^ eigenvalues u, two by two with opposite signs because 
of the AF symmetry. The solution for the GF Qn^n is given by 

IML 



9n,T 



IMI 



(22) 



with |M|,^ is the determinant made by replacing the n-th column of |M| by u in (1151) . 
Writing now 

\M\=Y[{uj-u,{Ky)), (23) 

i 

one sees that Ui (kxy) , ^ = 1, ■ ' ' ) 4A^2, are poles of the GF gn,n- Now, we can express gn,n 



as 



where fn {uji {Ky)) is 



Next, using the spectral theorem which relates the correlation function {S^Sj~) to the 
GF,— one has 

{S^S+) = lim^J J dKyJ_^ ^{9n,n' {uJ + ie) 

- 9n,n' {u^ - re)) ■ J^^^e^'^--^^^-''^ (26) 

where e is an infinitesimal positive constant and (3 = l/Zc^T, ks being the Boltzmann 
constant. 

Using the GF presented above, we can calculate self-consistently various physical quanti- 
ties as functions of temperature T. Large values of Ising-like interaction / will enhance the 
ordering. On the contrary, for / the transition temperature will go to zero according 
to the Mermin- Wagner theorem.— This is seen in the following. For numerical integretaion, 
we will use 80^ points in the first Brillouin zone. 

Figure [13] shows the sublattice magnetizations of the first two layers for = 2 with 
/ = —0.25 and / = —0.01 (upper and lower figures, respectively)). As seen, the surface 
sublattice magnetization is larger than the sublattice magnetization of the second layer for 

= 2 in qualitative agreement with the MC results shown in Fig. [71 in spite of the fact 
that due to a finite-size effect, there is a queue of the sublattice magnetization above the 
transition temperatures for MC results. Note that the AF coupling gives rise to a zero-point 
spin contraction at T = which is different for the surface spins and the second-layer spins. 

We show in Fig. [14] the phase diagram in the space (/, T) where Ti and T3 are transition 
temperatures of the surface sublattice 1 and the sublattice 3 of the second layer. 
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FIG. 13: Magnetization of sublattices 1 (surface) and 3 (second layer) versus temperature for 
= 2 and 7 = -0.25 (upper) and I = -0.01 (lower). 




FIG. 14: Phase diagram obtained by the GF method. Ti and are transition temperatures of the 
sublattice 1 of the surface and of the sublattice 3 of the second layer. The transition temperatures 
should go to zero as I ^ (see text). 

V. CONCLUDING REMARKS 



We have shown in this paper the crossover of the phase transition from first to second 
order in the frustrated Ising FCC AF film. This crossover occurs when the film thickness 



A^^ = 2 is smaller than a value between 2 and 4 FCC lattice cells. These results are obtained 
with the highly performing Wang- Landau flat histogram technique which allows to determine 
the first-order transition with efficiency. 

For Nz = 2, we found that in a range of temperature the surface spins stay ordered while 
interior spins are disordered. We interpret this as an effect of the frustration reduction: due 
to the lack of neighbors, the surface spins are less frustrated than the interior spins. As a 
consequence, interior spins are disordered at a lower temperature. This has been verified by 
the Green's function calculation. 

The second-order transition for A^^ = 2 is governed by the surface disordering and is 
characterized by critical exponents whose values are deviated from those of the 2D Ising 
universality class. We believe that this deviation results from the effect of the disordered 
interior spins which act as "correlated" random fields on the surface spins. We do not know if 
the critical exponents found here belong to a new universality class or they are just "effective 
critical exponents" which one could scale in some way or another to bring into the 2D Ising 
universality class. Anyway, these exponents seem to obey a weak universality. An answer 
to this question is desirable. 

One of us (VTN) thanks the University of Cergy-Pontoise for a financial support and 
hospitality during the final stage of this work. 
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